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ABSTRACT 

A finite element computer program has been developed for evaluating the 
transient behavior of surface temperature In high performance spur gears. The 
program Is an extension of a previously developed code that was limited to 
steady state behavior. The time dimension Is Implemented using two and three 
point finite difference schemes. The different schemes are provided to the 
user for the purpose of numerical stability and convergence studies. A de- 
tailed explanation of the gear cooling process leading to the establishment of 
a modified Blok model Is also Included. Other conventional models for approxi- 
mating the heat transfer coefficients are available for comparison. Prelimi- 
nary results are given showing snap shots of gear temperature contours at the 
Initial stages of tooth engagement. 

INTRODUCTION 

Amongst the three modes of gear failure, scoring Is known to be the most 
troublesome In high speed case hardened gearing. Experimental evidence [1] has 
shown that a sudden surge of the flash temperature of the gear takes place at 
the onset of scoring. Because of the Interdependence of the flash and bulk 
temperatures of the gear, an accurate prediction of the latter Is highly desir- 
able. Recently In Refs. 1 to 3 a finite element approach has been successfully 



applied to the estimation of bulk temperatures of a pair of meshing gears at 
steady state operation. The temperature distribution within each tooth Is 
assumed to be Identical and all transient effects are neglected. A variety of 
models were used to approximate the convective heat transfer coefficients with- 
in the elastohydrodyanmlc film attached to the tooth profile. The effects of 
oil Jet cooling and oil Jet Impingement depth were Included In Refs. 1 and 3. 
Because of the highly transient nature of thermal behavior at the onset of 
scoring a transient thermal analysis of the meshing gears would be extremely 
useful. Furthermore since the heat transfer coefficients, as will be shown In 
the text, are themselves time dependent functions, a steady state approach 
could not be used without major averaging approximations. 

The Intent of this paper Is to describe a computer analysis procedure for 
the transient thermal analysis of spur gears using the finite element method. 
The finite element equations are derived using the well known Galerkln or 
weighted residual approach as explained In detail In Ref. 4. The heat flux 
Input generated due to friction between the meshing gears Is derived as a time 
dependent function unlike the heat flux used In the steady state solutions of 
Refs. 1 to 3. Two separate models have been programmed Into the computer code 
for approximating the convective heat transfer process within the cooling film. 
The first allows the use of fixed heat transfer coefficients along the various 
sections of the gear profile as In Refs. 1 to 3 for the purpose of comparison 
with the steady state solution obtained therein. The second uses an array of 
time dependent heat transfer coefficients based on an asymptotic solution to 
the one-dimensional heat transfer equation as shown In the text. 
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Derivation of Finite Element Equations for Transient Heat Transfer In Gearing 
It Is shown In Ref. 3 that the majority of heat dissipation occurs at the 
tooth boundary. Therefore at the expense of little loss In accuracy, It can 
be assumed that the temperature within one tooth Is Independent of those within 
the neighboring teeth. This justifies the use of the one tooth sector model 
shown In Fig. 1, Instead of the whole gear. 

The partial differential equation governing heat transfer within the Iden- 
tified region Is given by 



T Current gear temperature, °C (°F) 

T Initial gear temperature, °C (°F) 

a 

2 2 

a Gear dlffuslvlty, m /sec (In. /sec) 
c Specific heat of gear material, J/kg °K (Btu/lb °F) 
k Thermal conductivity of gear material, W/m °K (Btu/ft °F) 
R q Outside radius of gear, m (ft) 

F Face width of gear, m (ft) 

q flv Average heat flux through the gear profile, W/hr (Btu/hr) 
The associated boundary conditions are 
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where 


n normal coordinate to the tooth profile 

2 

h^ convective heat transfer coefficient In tooth profile, W/hr-m °K 
( Btu/hr ft 2 °F) 

h $ convective heat transfer coefficient on front and back faces of gear 
P d diametral pitch of gear, m -1 (In. -1 ) 
q* q/q 

tooth Involute boundary, loaded side 
r. tooth Involute boundary, coast side 

t N 

r front and back end faces of tooth 
s 

Discretization In the space domain Is accomplished through the use of the 
Galerkln approach [4] In conjunction with the trial function. 

e(x, y, z, t) = l N^x, y, z) e^t) = N e (4) 

where are the spatial shape functions chosen as [2] 

N 1 = 2A (a 1 + b 1 +c i y)(1 " Bz2) (5) 

The Galerkln principle yields the following first order equation 


C6 + H6 + P = 0 


(6) 
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Proceedings In the manner described by Zlenklewlcz [4] for discretization In 
the time domain such that 


e(t) « e(t) = l L^t) e 1 

where L^(t) are time dependent shape functions and ^ represents at 
node n the time dimension. These are given as 


( 8 ) 
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0 < l < 1 


(9) 


where 


l = t/4t 

Applying the Galerkln principle to Eq. (6) by Integrating Its weighted average 
over a time element, It follows that 
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where Is any appropriately chosen weighting function of time (c). 

Noting the definitions of L n and L n+ ^ In Eq. (9), the following recur- 
rence relation for Is obtained. 

x n 

(C/At + H Y) e n+1 = [C/At - H(1 - Y>] ^ ♦ P (11) 

where 



( 12 ) 


If a shape function expansion similar to Eq. (8) Is used to express the heat 
Input term P, the last of Eq. (12) becomes 


! ■ £„ tl r * £ n (i - r) 

A variety of numerical schemes such as the Euler, Crank-NIcholson, backward 
difference etc. are the results of Eq. (11) for various assumptions of the 
weighting function Wj (see [4] for details). It can be shown that 0 < y < 1 
and that for numerical stability y > 1/2. 

Derivation of Time Dependent Frictional Heat Input 
The steady-state heat flux Input generated In friction between the two 
meshing gears was derived In [1] as: 



where the subscript 1 refers to the Instantaneous point of contact along the 
tooth profile. Because the analysis In [1] was restricted to steady state 
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running, the heat flux along the entire tooth profile was simultaneously ap- 
plied. In the present transient analysis, the Instantaneous flux q* at a 
time Instant t^ Is therefore given as 

q* = q^ i(t - t,) (14) 


where 6^ Is the Dirac delta function and q^ Is given by Eq. (13). 

The computer program Is designed to calculate, If and when the two gears 
are In contact, the corresponding roll angle measured from the beginning of 
engagement at the outer tip of gear and the corresponding radius of curvature 
p 1 . No heat Input Is generated beyond the lowest point of contact until the 
next time the tooth under consideration Is engaged. 

Convective Heat Transfer Models 

At the present time two basic options are available In the program. The 
first enables the user to use constant heat transfer coefficients Individual- 
ized In the manner described In [1]. These are presumably derived experimen- 
tally and averaged out over the expected temperature range of the gear tooth. 
In this model no convective heat dissipation Is assumed to take place beyond 
the fling-off angle and until the tooth Is reengaged. 

The other model for heat convection Is based on an asymptotic solution of 
the heat transfer equation along the tooth profile such that 


aj 

at 



(15) 


The exact solution of Eq. (15) for a constant temperature along the tooth pro- 
file gives the heat transfer coefficient at any Instant t as 


where K Is the thermal conductivity of the oil. Again In this model no heat 
dissipation Is assumed to take place beyond the fling off angle and until the 
tooth Is reengaged. 

The amount of heat flux dissipated In convection Is given by 

q Q = hA AT (17) 

where h may be assigned fixed values, while the tooth Is being cooled, or the 
program Is Instructed to calculate the current value of h from Eq. (16). 

RESULTS AND DISCUSSION 

The recurrence formula of Eq. (11) above constitutes a set of nonhomoge- 
neous algebraic equations for the dimensionless temperature array © n+1 at a 
time node n+1 In terms of the temperature array 6 n at the preceding time 
node n and the dimensionless heat flux array P. It should be noted that the 
latter Is made of a linear combination of the heat fluxes at times n and 
n+1, as defined In Eq. (13). Based on recommendations of well known finite 
element programs such as NASTRAN [5] the weighting factor y, convectlonal ly 
known as a finite difference parameter, Is assigned a value of 0.55. This 
value has been found to yield rapid convergence and greater accuracy for prob- 
lems with known exact solutions. The computational technique consists of a 
standard Initial value problem algorithm In which the Initial temperatures are 
assumed and heat flux vector Is computed from Eqs. (14) and (15). A Cholesky 
decomposition Is then administered to Eq. (11) to facilitate the solution [6]. 

The finite element grid used for computation Is shown In Fig. 2 In which 
the tooth Is divided Into 18 layers of triangular elements using the gear fi- 
nite element generator program ATOGEN. It was shown In Ref. 3 that In the 
steady state case, this mesh yielded results of sufficient accuracy. The 
example run Is that of the NASA four square testing apparatus described In 
Ref. 1. 
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The temperature distribution within the tooth Is conveniently represented 
In the computer output In the form of constant contour lines as Illustrated In 
Figs. 3 to 7. Each figure represents a chronological set of "snap shots" for 
the pinion tooth as It traverses the line of action between the highest and 
lowest points of contact. 

All the computer runs were made using an 8 pitch 8.89 cm (3.5 In.) pitch 
diameter, 0.635 cm (0.25 In.) wide spur gear operating at 10 000 rpm with a 
normal load of 5912 N/cm (3378 lb/in.) and a dimensionless oil Jet Impingement 
depth of 0.336. 

In Figs. 3(a) to (d) (Case 1) the heat transfer coefficients for the dif- 
ferent regions of the gear tooth model (fig. 1) are assumed constant during the 
period of contact and their values are Identical with those used In Refs. 1 and 
3 for the steady state temperature survey of the subject gear and were H $ = 0.6 
for the sides of the tooth, = 3.0 for the tooth flank surfaces without oil 
jet cooling, = 13 for tooth flank surface with oil Jet cooling and = 3.0 
for top land of tooth with oil jet cooling. As noted In Ref. 1, those values 
were obtained In an Iterative manner to yield best correlation with test data. 

The temperature contours are seen to be crowded around the point of con- 
tact, a would be expected, since the heat source at that Instant of time coin- 
cides with that point. The maximum temperature at any Instant during the mesh 
cycle always occurs at the Instantaneous point of contact. The contours pro- 
ceed to propagate and spread out down the tooth profile. It Is Interesting to 
note that closed loop contours are occasionally observed within the gear tooth 
Indicating the presence of local maxima In the thermal topology of the tooth.. 

' i 

These contours continue to propagate away from the contact point towards the 
unloaded coast side until they totally disappear. 
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The temperature contours at the same Instants of time are reproduced In 
Figs. 4(a) to (d) (Case 2) using a heat transfer coefficient hj obtained from 
the time dependent Eq. (16). Gear cooling Is constrained to the loaded side of 
the tooth above the lowest point of contact. The Initial distribution of tem- 
perature Is significantly different from the preceding case for the recess por- 
tion of contact, up to the pitch point, where the distribution approaches the 
case of steady heat transfer coefficients. In fact the maximum temperature 
rise during the contact period Is nearly the same for the two cases (14.18° C 
for case 1 versus 14.11° C for case 2). 

Case 3 (figs. 5(a) to (d)) plots pertain to cooling of the top land In 
addition to the loaded portion of the tooth. Comparison of Fig. 5 with Fig. 4 
Indicates no significant change with the additional cooling from the top land. 
One may therefore suggest that cooling the top land Is of little or no benefit 
In reducing the bulk temperature of the gear. More study Is required before 
this will be conclusive. The maximum temperature rise during the first engage- 
ment cycle for this case Is found to be 14.03° C. 

The remaining two cases consider the addition of cooling to the coast side 
of the tooth as well as the transverse planes of the gear. Because the tran- 
sient heat transfer formula, Eq. (16) Is Independent of temperature (except 
Implicitly through the thermal conductivity k and dlffuslvlty a), the In- 
stantaneous heat transfer coefficients used for all regions of the tooth are 
Identical. Although this Is far removed from practice. It Is Interesting to 
note the significant drop In temperature as a result of cooling the ends and 
coast sides of the tooth. The maximum temperature rise Is only a fraction of 
that computed for cases 1 to 3. The temperature contours are now continuously 
crowded around the contact point with most of the tooth remaining at the Ini- 
tial temperature. 
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An attempt was made to modify the heat transfer coefficients by multiply- 
ing them by a correction factor proportional to the steady state heat transfer 
coefficients, l.e., H $ x.6/13, H t x 3/13, Hj x 13/13 and H Jt x 3/13. The re- 
sulting contour plots are not shown In this paper because of their similarity 
to Figs. 6 and 7. 

The reduction In maximum temperature Is still very significant, although 
not as drastic as shown In cases 4 and 5 (figs. 6 and 7) (AT max = 6.03° C not 
shown versus 2.91° C (fig. 6(d)). 

Figure 8 Illustrates how the Instantaneous maximum temperature varies 
with time (or pinion rotation angle) such as would be detected In an experi- 
mental measurement. Although no experimental results are available at this 
time for the Initial stages of engagement, the wavy patterns obtained In 
Fig. 8 resemble those for the experimental data shown In Ref. 1. The minimum 
temperature noted In the graph Is clearly associated with the pitch point, 
where pure rolling takes place and friction losses are negligible. 

SUMMARY OF RESULTS 

A finite element computer program has been developed for evaluating the 
transient behavior of surface temperature In high performance spur gears. The 
program Is an extension of a previously developed code that was limited to 
steady state behavior. The time dimension Is Implemented using two and three 
point finite difference schemes. The different schemes are provided to the 
user for the purpose of numerical stability and convergence studies. A de- 
tailed explanation of the gear cooling process leading to the establishment of 
a modified Blok model Is also Included. The following results were obtained: 

(1) The maximum temperature obtained with constant heat transfer coeffi- 
cients, used to predict steady state conditions, Is nearly the same as that 
obtained from a time-dependent modified Blok model when heat dissipation Is 
constrained to the loaded side of the tooth. 
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(2) Cooling of the top land was found to be of no benefit In reducing 
gear temperature. 

(3) Cooling of coast side of the tooth was, on the other hand, of great 
significance and should be utilized whenever possible. 

(4) The temperature rise along the recess portion of the gear tooth Is 
markedly lower than the approach portion. It should be advantageous , there- 
fore, to use recess-action type gears to preclude scoring. 
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(a) Pinion rotation angle, 0.72°; time, 0. 012xl(f 3 seconds from 
point of initial tooth contact 



(c) Pinion rotation angle. 13.68°; time, 0. 228 xl(f 3 seconds from 
point of initial tooth contact. 




(b) Pinion rotation angle, 5.04°; time, 0. 08x 10“ 3 seconds from 
point of initial tooth contact 


D O 

a o 

(d) Pinion rotation angle, 2a 88°; time, a 348x10^ seconds from 
point of initial tooth contact 


Figure 3. - Transient temperature isotherms using steady state heat transfer coefficients H s , a6; H t and H: t , 3; H:. 13. For 8-pitch spur gear: pitch diam- 
eter, 8189 cm (3. 5 in. ); face width, a 635 cm (a 25 in. ); speed, lOOOOrpm; load, 5912 N /cm (3378 lb/ in. V, dimensionless impingement depth, a 336. 
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(a) Pinion rotation angle, 0.72°: time, 0.012xl(T 3 seconds from 
point of initial tooth contact 
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(cl Pinion rotation angle, 13.68°; time, 0. 228 xl0“ 3 seconds from 
point of initial tooth contact. 



point of initial tooth contact. 



Figure 4. - Transient temperature isotherms using transient thermal analysis with computed heat transfer coefficients, load side cooling only. For 8-pitch 
spur gear: pitch diameter, 8. 89 cm (3.5 in.); face width, 0.635 cm (0. 25 in.); speed, 10 000 rpm; load, 5912 N/cm (3378 Ib/in. ); dimensionless impinge- 
ment depth, a 336. 
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(a) Pinion rotation angle, (X 72°; time, Q.012xl0“ 3 seconds from 
point of initial tooth contact. 
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(c) Pinion rotation angle, 13.689; time, Q,228xl0~ 3 seconds from 
point of initial tooth contact. 


(b) Pinion rotation angle, 5.04; time, 01 08 x 10" 3 seconds from 
point of initial tooth contact 




Figure 5. - Transient temperature isotherms using transient analysis with computed heat transfer coefficients, load side and top land cooling only. 
For 8-pitch spur gear: pitch diameter, & 89 cm (3. 5 in. ); face width, G. 635 cm (a 25 in. ); speed, 10 000 rpm; load, 5912 N/cm (3378 Ib/in. ); 
dimensionless impingement depth, a 336. 
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(a) Pinion rotation angle, 0,72°; time, a 012xl0' 3 seconds from (c) Pinion rotation angle, 13.68°; time, a 228xl0~ 3 seconds from 

point of initial tooth contact. point of initial tooth contact. 
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(b) Pinion rotation angle, 5. 04°; time, a 08xl0" 3 seconds from (d) Pinion rotation angle, 2a 88°; time, a 348xl0~ 3 seconds from 

point of initial tooth contact point of initial tooth contact 

Figure 6. - Transient temperature isotherms using transient analysis with computed heat transfer coefficients modified by the ratio of the coefficients 
used in (ref. 1). Cooling of the loaded side, coast side and top land only. For 8-pitch spur gear: pitch diameter, & 89 cm (3. 5 in. ) ■ face width, 
a 635 cm (0.25 in.); speed, 10 000 rpm; normal load, 5912N/cm (3378 Ib/in. ); dimensionless impingement depth, 0.336. 
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(a) Pinion rotation angle, a 72°; time, a 012 xlO" 3 seconds from 
point of initial tooth contact 


(c) Pinion rotation angle, 13.68°; time, 0,228xl0' 3 seconds from 
point of initial tooth contact 




(b) Pinion rotation angle, 5.04°; time, a08xl0 -3 seconds from 
point of initial tooth contact 


(d) Pinion rotation angle, 2a 88°; time, a348xl0‘ 3 seconds from 
point of initial tooth contact 


Figure 7. - Transient temperature isotherms using transient analysis with calculated heat transfer coefficients modified by the ratio of the coefficients 
used in (ref. 1). Cooling of the loaded side, coast side, top land and end of the teeth. For 8-pitch spur gear: pitch diameter, a 89 cm (3.5 in. ); 
face width, a635 cm (a 25 in. ); speed, lOOOOrpm; normal load, 5912 N/cm (3378 lb/in.); dimensionless impingement depth, a 336. 
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